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Abstract 

Recent experimental studies of ^''Rb spinor Bose Einstein condensates have shown the existence 
of a magnetic field driven quantum phase transition accompanied by structure formation on the 
ferromagnetic side of the transition. In this theoretical study we examine the dynamics of the 
unstable modes following the transition taking into account the effects of the trap, non-linearities, 
finite temperature and dipole-dipole interactions. Starting from an initial state which includes 
quantum fluctuations, we attempt to make quantitative comparisons with recent experimental data 
on this system and estimate the contribution of quantum zero-point fluctuations to the domain 
formation. Finally, using the strong anisotropy of the trap, we propose ways to observe directly 
the effects of dipole-dipole interactions on the spinor condensate dynamics. 

PACS numbers: 03.75.Mn 
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I. INTRODUCTION 



The physics of phase transitions between ordered and disordered phases has been studied 
extensively in the past and has yielded answers to many fundamental questions about the 
effects of interactions and of quantum and thermal fluctuations on the equilibrium properties 
in various phases. However, because of the short time-scales and high levels of noise and 
impurities involved in traditional condensed matter systems, it has been difficult to study 
the dynamics of such transitions in detail. Ultra-cold atoms, because of their low densities 
and temperatures, provide an opportunity to study the non-equilibrium dynamics around 
phase transitions in a spatially and temporally resolved fashion. Bose Einstein condensates 
(BECs) of atoms with a spin degree of freedom, i.e. spinor BECs, are an example of such 
a system where it is possible to observe non-equilibrium spin dynamics and how they are 
affected by quantum noise and the proximity to phase transitions^. 

Recent experiments^-'^ have reported the formation of magnetic structures in ultra-cold 
spin-1 *^Rb gases. Spin-1 atoms support two characteristic families of quantum states: polar 
states exemplified by the = hyperfine state and magnetic states exemplified by the 
iriz = ±1 hyperfine state^'^'^; here denotes the eigenvalue for the dimensionless spin 
projection along the z axis. As discussed further below, spin-dependent contact interactions 
naturally favor the ferromagnetic state in *^Rb spinor condensates. The Rb atom may also 
be subjected to an extrinsic quadratic Zeeman shift, which lowers the energy of the = 
state by an amount q with respect to the average of the energies of the = ±1 states. 
In the experiments, a non-magnetic = condensate is prepared at a large value of q, 
where its internal state composition is stable. Following this, q is rapidly quenched to a 
regime where the initially prepared = becomes unstable. As population flows to the 
initially unoccupied = ±1 states, the condensed gas is observed to break translational 
and rotational symmetry spontaneously and form ferromagnetic domains of transversely 
magnetized atoms. In the experiment, the system also evolves under a signiflcant linear 
Zeeman splitting of the three Zeeman sublevels = 0, ±1. However, as we shall see later, 
since most of the terms of the Hamiltonian are invariant under global unitary spin rotations, 
the linear Zeeman shifts can be eliminated from the theoretical treatment by a global unitary 
transformation of the spin into the rotating frame of the Larmor precessing spin. 

The symmetry-breaking ferromagnetic domain formation discussed above is a phe- 
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nomenon that accompanies a large number of thermodynamic phase transitions such as 
the paramagnetic to ferromagnetic phase transition seen in iron when the temperature is 
lowered below its Curie point. However, unlike the thermodynamic phase transition in Fe, 
the symmetry-breaking transition in the spinor BEC is a quantum transition that can occur 
at arbitrarily low temperature. Concomitantly, the initial fluctuations that seed the sym- 
metry breaking dynamics in a spinor BEC need not be thermal in nature but, rather, may 
be of a quantum origin. The quantum origin of such fluctuations makes the detailed study 
of the origin and dynamics of the spontaneous magnetization of fundamental interest. 

One approach to developing a theoretical understanding of the above phenomenon is to 
analyze the low energy dynamics of a spinor BEC by linearizing the Heisenberg equation of 
motion of the bosonic annihilation operators for the atoms. In previous theoretical studies^i^, 
these linearized equations of motion were obtained by expanding the fields corresponding 
to the three components of the BEC around the initial = state. One finds that the 
low energy dynamics of the condensate in the initial state can be described by a set of three 
low energy excitations composed of a gapless phonon mode and two gapped magnon modes. 
On quenching the quadratic Zeeman shift to a value which is below the phase transition, 
the magnon modes are found to become unstable. These unstable modes amplify any initial 
perturbations from a homogeneous polar state. Thus the domain formation following the 
quench is described as resulting from the quantum fluctuations in the initial ground state 
being dynamically amplifiecl^'^. 

The theoretical calculations discussed above treat the domain formation by calculating 
the linearized dynamics of a homogeneous spinor condensate with local interactions. While 
these calculations yield results in qualitative agreement with experiment, a quantitative com- 
parison to experiment is essential for confirming the quantum nature of the initial seed and 
of the amplifier driving the structure formation. An improved understanding of the dynam- 
ics of the spinor condensate including the effect of the trapping potential and non-linearities 
has been obtained in previous theoretical worka^"^*^. In the current work, we improve fur- 
ther our understanding of the dynamics of spinor condensates by including the effects of 
dipole-dipole interactions and finite temperature and use this improved understanding to 
determine the domain formation resulting from quantum fluctuations within the truncated 
Wigner approximation (TWA). 

We begin with a discussion of how the condensate in a pancake-shaped trap can be 
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modeled by a two dimensional Hamiltonian. Next we introduce a general framework for 
calculating the eigenmodes of the inhomogeneous gas which are used to describe the time 
evolution of the initial quantum fluctuations. Following this we discuss explicitly the effect 
of dipole-dipole interactions, the trap and finite temperatures on the dynamics of a spinor 
BEC. Next we introduce the approximations and the general computational framework that 
allows us to include all these effects together with non-linearity effects such as saturation of 
the transverse magnetization at long times. Finally we comment on how our results compare 
to experimental results. 



II. TWO-DIMENSIONAL EFFECTIVE HAMILTONIAN 

The many-body Hamiltonian for the spin-1 ^''Rb gas considered in this work, expressed 
using bosonic fields -ipTny to represent the three magnetic sublevels irip = 0, ±1, is given asr^ 

H = J dR 1^ 5Z 1^^"!' + ~ ^) + ^2F' + qit)Q + Ft,ap(R)n| + Udipoie (1) 

where M is the atomic mass, and the quadratic Zeeman shift, q{t) is applied along P 
and Vtrap(R) is the external trap potential. The parameter cq is the strength of the spin- 
independent atom-atom short range repulsion, C2 is the coupling constant for the spin- 
dependent contact interaction, and Udipoie is the dipole-dipole interaction. The boson field 
operators appear in the magnetic part of the interaction terms through the spin density 
operators defined as = 7 ^a,/3,7'^/3'^7 where Ja are the spin-1 matrices in the fun- 
damental representation in a basis such that the spin is quantized along P, which is the 
axis determined by the quadratic Zeeman shift. The quadratic Zeeman shift term in the 
Hamiltonian is of the form Q = J^/s-yi^ ' JYp^-y'^^'^'y The total atom density is given by 
h = I'^/'aP- In the above discussion, the condensate is prepared in the rrip = hyperfine 
state, where rup is the magnetic quantum number for spins quantized along the axis P. This 
state is prepared at a high quadratic Zeeman shift q{t) ^ 2|c2|n3/) where 1130 is the peak 
density at the center of the condensate. The quadratic Zeeman shift is then rapidly reduced 
below the critical value of 2|c2|n3£) and the state of the condensate is allowed to evolve. In 
the experimental work discussed in this paper the axis along which the quadratic Zeeman 
shift is applied P is taken to be the z-axis. 

As mentioned in the introduction, the atoms in the system evolve under a strong magnetic 
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field which we have eliminated from the Hamiltonian described in Eq. 1 , by transforming to a 
rotating frame. Such a transformation leaves the rotation-invariant terms in the Hamiltonian 
unchanged but affects the quadratic Zeeman shift term and the dipole-dipole interaction term 
Udipoie- Yet, under the experimental conditions that the Larmor precession frequency is far 
higher than that of the interaction-driven spin dynamics, one may consider the quadratic 
Zeeman shift term and the dipole interaction as precession-averaged static terms in the 
Hamiltonian. It may also be possible to vary dynamically the axis of the quadratic Zeeman 
shift to follow the Larmor precession of the initial state applied by the initial magnetic 
field, in which case the quadratic Zeeman shift term would be intrinsically stationary in the 
rotating frame of the spin. 

The above defined Hamiltonian can be used to calculate the dynamics of the transverse 
magnetization, which is the observed quantity that is used in the above mentioned experi- 
mental studies to observe the non-magnetic to ferromagnetic phase transition^. The x and 
the y components of the transverse magnetization, and Fy, can be combined into a single 
complex transverse magnetization operator Fj_(R) = F^^R) +Fy{n) which is given in terms 
of the fields ipa by 

F^(R) = V2 (^5(R)V;+i(R) + (2) 

The random magnetization domain pattern that forms after the quench can be characterized 
by a correlation function of the above defined transverse magnetization, which we take as 
^(^R) = / dR ^/1(R-(5R/2)F|(R + (5R/2)^/ (/ dR(n(R)))^l. 

An important feature of the experiments is the use of condensates with widths in one 
dimension (y) that are smaller than the spin healing length, 27r/ ^y2M\c2\nsD- The three- 
dimensional Hamiltonian is thus reduced to a two dimensional form through the substitutions 
■^a(R) = 0Q,(r)^(r; y)^n2D where ?/ = R ■ y and R = r + yy. The peak density integrated 
along the y axis is given by n2D = / dy (^^(O, y))initiav Here ^(r; y) represents the normahzed 
spatial profile of the wavefunction at each point r in the two dimensional x-z plane. The scale 
of the initial quantum zero-point fluctuations of the two dimensional field, 0^; is determined 
from the canonical commutation relations — 'i^^^^ ~ ^')^a,i3- 

Using these relations and the Heisenberg equations of motion one can construct a time- 
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evolution equation for the operator 0a (r) of the form 

'd,k{v) = V^4W +4W j ds {-^ {Vi(r;y))i{v;y) 

+ ^„(r)/d!,«^(r;,)V;„,(r; + (3) 

The second term on the right hand side of the above equation may be considered to be an 
effective renormahzation of the potential energy related to confinement effects. Within the 
Thomas- Fermi (TF) approximation this term is small through most of the condensate and 
is hence ignored. 

For simplicity, we assume the wavefunction profile to be a TF profile given by ^(r;y) = 

^5^(1 - yVi?|p,^^(r)) V2 where i?TF,,(r) = i?TF,,(l - x^lR\^^., - ^V^^If, J'^' and i?TF,x, 
i?TF,j/ and /2tf,^ are the TF radii in the x, the y and the z directions respectively. With the 
above choice of a transverse profile, the contribution of a local two-body potential of the 
form c^5(r — r')S{y — y') to the interaction term Uinteraction simplifies to 

Cm f ^'m{r;y)dy = c^T^, I e{r;y)dy<p'{r) = Q.8c^n,n^^<P'{r) (4) 

where the index m is either or 2 depending on whether we are referring to the spin 
independent or spin dependent parts of the contact interaction, respectively. In the rest of 
the article we will be using a two-dimensional position-dependent effective interaction by 

Cm(r) = 0.8c^n3i3i?TF,j/(O)/-RTF,2/(r). 



III. QUANTUM DYNAMICS AND QUANTUM NOISE SEEDED DOMAIN FOR- 
MATION 

Let us now describe the fluctuations and domain formation in terms of the two dimen- 
sional fields derived above. We consider the dynamics and low energy fluctuations of the 
initial nip = state by shifting the operator corresponding to the nip — component by 
00 — + Vo, where n(r) = (0o(r)0o(i")) is the equilibrium density for the condensate 

at large positive quadratic Zeeman shift of ^ ^ |c2|7T'3I3. The Hamiltonian can now be ex- 
panded to second order in the small fluctuations in the small fluctuation operators 0±i,57o- 
In this Hamiltonian, the terms involving 170, which describe the scalar Bogoliubov spectrum 
of phonons and free particles, separate from those involving spin excitations; these latter 
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terms provide the following Hamiltonian 

n I 



Hmagnon = / 0+lV 0+1 + 0LiV 0-1 



+ J {q{t) + C2(r)n(r) + Vtrap{r) + Co(r)n(r)) (^^i^^+i + '^-i^-i) 

+ [ C2(r)n(r) + ^^i^li) • (5) 




For simplicity, the effect of the dipole-dipole interaction term has been ignored here and 
its discussion is postponed to Section V. The dynamics of the magnetic degrees of freedom 
obtained from the above approximate Hamiltonian are given as 

dt^ir) = (^-^ V + + Vtrapir) + (co(r) + C2(r))n(r)^ a, + C2n(r)a^| <l>(r). (6) 

where we have introduced the spinor $(r) = 

The above spinor equation of motion can be used to describe the dynamics of the con- 
densate in terms of normal modes T^^\r) with frequencies ±i?„. The dynamics of $(r, t) 
are then determined as 

n,CT=itl 

where d^T^ are the mode occupancy operators. The magnon modes are stable when the 
eigenenergies En are real, and unstable when En are complex. It is these unstable modes 
that amplify quantum fluctuations to generate macroscopic magnetization in the quenched 
spinor gas. 

In the case of a homogeneous condensate the normal modes Tl'''^''(r) can be reduced to 
a product of a plane wave state and the momentum dependent two-component spinor that 
appears in standard treatments of the linearized Gross-Pitaevskii equations. However the 
determination of these modes in the case of an inhomogeneous density must be done in real 
space using explicit numerical diagonalization of a generalized eigenvalue problem. In the 
case of a positive quadratic shift, which is the focus of in this article, the frequencies of these 
eigenmodes can be shown to be either purely real or imaginary as discussed in Appendix 
A. A similar eignmode expansion for a trapped spinor condensate in the limit of vanishing 
quadratic Zeeman shift has been reported in previous work-. 

The quantum noise amplified by these unstable modes is entirely contained in the cor- 
relation function of the spinors relative to the initial state. Since the initial state of our 



system is assumed to be prepared as a condensate of atoms in the rup — state, the 
population in the nip — ±1 states is neghgible and the relevant correlator is given by 
^$(ri, 0)l>^(r2, 0)^ = 5(ri — r2)(l+(T2)/2. In the description of the dynamics of the conden- 
sate in terms of magnon modes, the quantum fluctuations become encoded in the quantum 
mode occupancy operators di^\ which can be derived from the spinor operator at the time 
of quench l>(r, 0), dn^ — J dr Ti'^''*(r)l>(r, 0), where Ti'^'*(r) are the dual modes, explicit 
expressions for which can be found in Appendix A. 

Given this initial noise and the linear dynamics of these magnon modes, we may calculate 
the magnetization correlations that may be observed at short times after the quench. Lin- 
earizing the transverse magnetization as F±{r) = n2D \/2n(r) (j)\_i{v) + 0_i(r)j, we obtain 
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_ / Mhjr - 5r/2)F[ir + Sv/2)) _ ^ ^.(.^„_.^^^), 
{J drn2Dn{r)Y 

/ dr^n{r - 5r/2)n{r + 5t/2) T^;:')^(r - 5v/2){l + a^)T^r^(r + 5r/2) 



(/ drn(r)) 

X j drit(f)t(rO(l + Ot(„-')(rO (8) 

Note that the correlation function defined above suffers from a UV divergence. The 
physical origin of this divergence is the fact that F±_ represents the magnetization of point- 
hke particles. This UV divergence is however not observed experimentally because of the 
physically natural cut-offs like the finite spatial and temporal resolution of the measuring 
apparatus, which introduces a natural spatio-temporal averaging into the observed magne- 
tization. For our purposes it suffices to consider only the contribution of unstable magnon 
modes to the transverse magnetization. Hence in our calculations we avoid this UV diver- 
gence in the correlator of F±_ by restricting the above sum to imaginary frequency modes. 



IV. NON-LINEARITY EFFECTS: TRUNCATED WIGNER APPROXIMATION 

In the last section we saw how we can describe the physics of the quench by expanding 
the Heisenberg equations of motion about the initial condensate state and keeping terms up 
to linear order in the fluctuations. Even though this might be expected to be a relatively 
accurate description at short times where deviations from the initial state are small, it breaks 
down at longer times and predicts an unphysical diverging magnetization. 
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Such a divergence is avoided by considering the complete Hamiltonian, including the 
higher order terms neglected in our prior approximation. The direct solution of the Heisen- 
berg equations of motion with the non-linearity terms would present an extremely difficult 
task. This is a general feature of problems involving quantum many-particle systems and for 
this reason various approximations must be used to understand such problems. One such 
approximation, that has been seen to describe the dynamics of BECs reasonably welUi, is 
the Truncated Wigner Approximation (TWA)^^. 

Within the TWA one is interested in calculating the time evolution of the Wigner Distri- 
bution Function (WDF) of the fields ^a{'^), which are the quantum analog of the classical 
phase space distribution. The fields are assumed to evolve according to the Gross-Pitaevskii 
equation (GPE) which is the mean field version of the Heisenberg equation of motion where 
the operator ^^(r) has been replaced by the time- dependent order parameter 0o-(i')- Quan- 
tum fluctuations around the mean field time evolution are included in the TWA by adding 
to the initial condition of the order parameter a noise term which is picked randomly from a 
classical distribution of fields corresponding to the initial WDF. In our case, where the initial 
state has a negligible population in the rrip = ±1 hyperfine state, the classical distribution 
of the randomly picked wavefunctions, 0o-(r), is a Gaussian distribution with variance given 
by {4>l{r)(j)^'{r')) classical = ^(^iratia«|{0l(i")0f7'(r')}|^miHai) whcrc {'^initial) is the initial state 
with all atoms in the rrip = hyperfine state.—. 

Within the TWA, we obtain the magnetization correlation function by averaging over 
numerical results obtained for the different, random representations of quantum noise. As 
shown in Fig. 1, this procedure yields a satisfactory result for the magnetization variance 
G{0) that saturates rather than diverging. 

The TWA has been shown to describe dynamical phenomena in BECs with reasonable 
success. However the TWA fails to describe certain aspects of the dynamics of BECs, such 
as the non-condensate fraction, and improvements beyond the TWA have been proposed in 
several works^*^. It is difficult to estimate the accuracy of the TWA in the full multi-mode 
spinor condensate system that we are studying. However as has been observed in previous 
work^'^, it is possible to solve the Hamiltonian under consideration within the single mode 
approximation exactly. This corresponds to the limit of a small trap where we can ignore the 
spatial dependence of the dynamics of the atoms completely and the system of 2N atoms can 
be described by a Hamiltonian given by Hsma = (2A^— 0j0o)(0+i0+i+0li0-i) — (0+i0-i0o+ 
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FIG. 1: Evolution of the variance of the transverse magnetization for the experimental 
configuration^ for q/h = 2 Hz and C2(0)//i = 8 Hz, calculated using the TWA to include non- 
linearity induced saturation effects. 

The time evolution of the fluctuations in the transverse magnetization of this 
Hamihonian can be determined exactly by numerically solving the time evolution of an initial 
state where all atoms are in the mp = hyperfine state. We compared the results of this 
calculation for a system of 2000 atoms to the time evolution of the transverse magnetization 
obtained within the TWA for the same system and found excellent agreement between the 
growth rate within the TWA to exact results up to the saturation time within the TWA. 
The exact transverse magnetization was found to increase past the saturation value within 
the TWA to a value that was 10% higher than the TWA saturation. While it is possible 
that the multimode non-linearity of our system causes physics beyond the TWA to become 
directly relevant, the above comparisons of the TWA to exact results for the single mode 
systems demonstrate that the TWA accounts for some of the effects of non-linearity in these 
systems. 

V. FINITE TEMPERATURE EFFECTS FROM AN INITIAL PHONON POPU- 
LATION 

In attempting to make quantitative comparisons between calculations and experimental 
observations, it is imperative to consider the role of the non-zero temperature on the initial 
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preparation and later evolution of the experimental system. In fact, at first glance, one 
might expect the quench experiments reported to be wholly dominated by thermal effects, 
given that the gas is prepared by evaporative cooling at a temperature of T ^ 50 nK for 
which the thermal energy is far larger than the spin-dependent energies responsible for the 
quench dynamics, i.e. fc^T ^ |c2|n3D- However, one must consider separately the kinetic and 
spin temperatures of the paramagnetic condensate in these experiments. While the thermal 
population of the scalar excitations, the nip = Bogoliubov excitations about the rup = 
condensate is indeed determined by the T = 50 nK kinetic temperature of the gas, the 
magnon excitations are expelled from the gas by the application of magnetic field gradients 
that purify the nip = atomic population. To the extent that such state purification is 
effective, and that magnon excitations are not thermally produced, e.g. by incoherent spin- 
exchange collisions, the initial spin temperature of the system is indeed near zero. Thus 
remarkably, a purely quantum evolution may indeed occur in the non-zero temperature gas. 

Here, we consider the possible influence of the thermal population of scalar excitations 
on the quantum quench experiments. The study of the coupling of phonons and magnons 
requires going beyond the linearized Heisenberg equation of motion. Thus, we consider 
the time evolution operator for the quantum state of the spinor condensate as a coherent 
state path integral U(ti,t2) = f Y[a^'PoiD(f)*^exp{zS[(j)a,4>a]), as has been found useful for 
many boson problems^^. Here S is the action for the three-component boson field corre- 
sponding to the Hamiltonian in Section II. The scalar phonon fluctuations are composed 
of a scalar density fluctuation, 6n{r,t) = J2a\'^o'i^j^)\^ ~ ''^(i")? and current fluctuations 
associated with the density fluctuations, required by number conservation. In the case 
of a condensate with population dominantly in the = hyperfine state, the current 
fluctuations can be described by the superfluid phase A(r, t) defined through 0o(r, t) = 
e*'^*^'"'*^ ^/n{r) + 5n{r, t) — |0+i(r, t)\'^ — |0-i(r, Since the density fluctuations are gapped 
at a high energy by the term co5n(r,t)^/2 in the action S, they can be integrated out to 
leave an effective action involving the phase A(r, t). Therefore in order to eliminate the high 
frequency density fluctuations we perform the field substitution 0±i(r,t) (f)±i(r,t)e'^^^^'*^ 
in the action S and then integrate out the density fluctuations 6n{r,t). This leads to 
the approximate action Sapprox = Sphonon[>^] + Smagnon[4>±i\ + Vcoupiing[(l)±i, M whcrc the 
Sphonon = nX^/2co + n^(VA)^, Smagnon is the usual action for the (j)±i atoms without 
the scalar interaction term and VcoupUng = ' Xla ■'■^(^a^^a)]- The last term is the 
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interaction that describes the couphng between phonons and magnons. 

The non-zero kinetic temperature of the gas gives causes (low frequency) fluctuations of 
superfluid phase with variance given by (|Afcp) = 2MkBT/h'^k'^n2D- These thermal phase 
fluctuations couple to the dynamics of the magnons, through the interaction term VcoupUng- 
A rough estimate of the magnitude of the effect of the kinetic temperature can be made by 
considering the dimensionless ratio of the r.m.s value of VcoupUng to the spin mixing energy 
\c2\n which is given by ^JWk^~j2M ^ {kBT)klf^^^^Jn2D/\c2\n^D ~ \/ {2MkBT) | C2 1 / h^n2DCQ 
where kspm is the wavevector associated with the spin healing length and kphonon is the small 
wavevector associated with a phonon at the energy scale of the spin dynamics, |c2|n3£). For 
the kinetic temperature in experiment of 50 nK, this dimensionless parameter characterizing 
thermal effects is found to be less than 1.3 x 10~^. A more rigorous calculation within the 
TWA, where phonons are introduced by adding random thermal fluctuations to the initial 
conditions in (pQ, conflrms our rough estimate by showing a negligible effect of the kinetic 
temperature. 

VI. ROLE OF DIPOLE-DIPOLE INTERACTIONS 

In the preceding paragraphs we have discussed the physics of the formation of domains 
from quantum fluctuations in a trapped quasi-two-dimensional condensate with ferromag- 
netic interactions. However theoreticay^*^ and experimental studies^ suggest that dipolar 
interactions play an important role in determining the magnetization textures for this sys- 
tem. In this section, we provide the flrst characterization of the role of dipolar interactions 
on the quantum quench dynamics of a ^^Rb spinor BEC. 

The atomic spin undergoes Larmor precession at a high frequency, on the order of tens of 
kHz, even as slower dynamics responsible for spontaneous magnetization transpire. While 
this Larmor precession has no influence on average on the spin dependent s-wave contact 
interaction or the quadratic Zeeman shift, the time averaged Larmor precession of the atoms 
must be accounted for in calculating the influence dipolar interactions, yielding an effective 
precession-averaged interaction of the form^i^ 

n ^2 f .J, (Ri - R^)^ - 3(D ■ (Ri - R^))^ 

Svr J IR1-R2P 

[3(D • F(Ri))(D ■ F(R2)) - F(Ri) ■ F(R2)] (9) 
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where D is the dipole-precession axis (the magnetic field axis), (yf^ = 1/2 is the gyromagnetic 
ratio of the electron, /x^ is the Bohr magneton. Integrating over the thin dimension of the 
condensate, we derive an effective two-dimensional dipole interaction as 

Cfipol. = ^ f dridr2A'(ri,r2)[3(D.F(rO)(D.F(r2))-F(ri)'F(rj)] (10) 



A-(r.r') = IR R.^. ' e(r;y)eir':v')dydv' (11) 

where the dipole interaction strength is given by Cdd = ^^^J^" (grfJ'B)'^- The dipole-dipole 
interaction term Udipoie is a spin-dependent interaction term in the Hamiltonian in addition 
to the ferromagnetic part of the contact interaction already discussed. The total of the two 
spin-dependent parts of the interaction Hamiltonian is given by 

Hsp^n= I dvC2{r)¥{vf + Ud^pole= [ dn dr2 cf (n, Ta) ^-^^""'^ ' ^^""'^^ 



A/ra(ri)n(r2) 

+ ^ j dridr2Qrfi^(ri,r2)(D-F(ri))(D-F(r2)) (12) 
cf (n, ra) = C2(ri)n(ri)5(ri - r^) - Cdd^M^iH^Kir^, r2)/2 (13) 

Thus apart from renormalizing the spin-dependent part of the contact interaction to d^, 
the dipole-dipole interaction also has an intrinsically anisotropic contribution which is given 
by the second term in Eq. 12, where the anisotropy is not related to the spatial anisotropy 
of the dipole interaction kernel K. This term however turns out to not be relevant for 
the linearized dynamics in the case where the dipole precession axis D coincides with the 
spin-quantization axis P. 

In the homogeneous case c"^ can be written as below 

cfik,x,v) = C2iO)-^Kik,x,v) (14) 

where K{k, x, v) is the dipole interaction kernel, and x is the angle that D makes with the 
y-axis and rj is the polar angle of the vector D in the plane of the BEC, as shown in Fig. 2. 
The wave-vector k is taken to be along the z-axis in the plane of the BEC. As discussed in 
Appendix A, the eigenmode treatment discussed in Section III can be easily generalized to 
include dipole-dipole interactions. 

From Eq. [14] it is apparent that in the three-dimensional homogeneous case dipole-dipole 
interactions enhance structure formation for wave-vectors along the dipole-precession axis 
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FIG. 2: Momentum dependence of the contribution of the dipole interaction kernel K defined in 
Eq. [15] to cl^ for Cdd/h = 0.8 Hz, C2{0)/h = —8.0 Hz. The inset shows the orientation of the 
dipole-precession axis D relative to the coordinates and the plane of the BEG. 

and suppress it for wave-vectors transverse to the dipole-precession axis. However the effect 
of dipole-dipole interactions on a quasi- two-dimensional condensate is qualitatively different. 
The Fourier transform of the interaction K{k, x, v) the case of the parabolic TF transverse 
profile along the y direction is difficult to compute analytically. To obtain a qualitative 
understanding of dipole interactions we consider the case of a 2D condensate for the case 
of a Gaussian profile C,{r;y) of width w = ^-^i?xF,j/- This width is chosen so that the 
peak density for the normalized profile matches that of the TF profile. The expression used 
for the Gaussian regularized dipole interaction derived in the Appendix B can be used in 
conjunction with standard integrals to determine the Fourier transform of K{r) for this 
Gaussian choice of profile to be 



^^(1 - -sin^x) + 7Tkw{sin'^xcos2r] + Ssin^x - 2)(1 - Erf(A;u;))e'=''" 
3 2 

(15) 

The contribution of the momentum variation of the dipole interaction kernel K to cf^ is 
shown in Fig. 2. In the large k limit, this expression, apart from a factor of \/2 arising from 
the re-normalization because of the transverse profile, reduces to —^(1 — 3 sin^ x cos^ 77) 
which is the three dimensional form as expected. As seen in Fig. 2, the small k limit is 
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FIG. 3: Magnon dispersion curves for an unbounded two-dimensional condensate in the x — z 
plane, for q/h = 2 Hz, C2(0)//i = 8 Hz and Cdd/h = 0.8 Hz, including the effects of dipole-dipole 
interaction for D aligned along x,y or z, and the wave-vector being assumed to be aligned along z. 

found to be consistent with previous theoretical studies^. 

In the case where the dipole-precession axis D coincides with the spin-quantization axis 
of the atoms P, we can use the explicit form for given in Eq. [T3] to discern the effect 
of the anisotropic dipole-dipole interactions on behavior of the spin dynamics by studying 
the dispersion relation -E^(k) = (||^ -|- q)i^ + q + 2c'if(k)) in the presence of dipole-dipole 
interactions. As seen in Fig. 3, when the dipole-precession axis points along the long axis of 
the condensate, i.e the 2;-axis, as in the experiments, the effect of the dipole interaction is 
weak and the dipole interactions slightly shorten the length scale and lengthen the time scale 
of domain formation. In contrast, the domain formation is significantly slowed down by the 
dipole-dipole interaction when D is oriented along the x-direction. Interestingly when the 
dipole axis is pointed along the y-axis, the thin axis of the condensate, the rate of structure 
formation is dramatically increased in both directions. 

The effects of the anisotropic dipolar interactions may also be highlighted in quantum 
quenches where the dipole-precession axis D differs from the spin-quantization axis P. 
Specifically, consider the case where D = while P is prepared to be orthogonal to D (i.e. 
the axis P Larmor precesses in the x — z plane). In this case the intrinsic spin-anisotropic 
term in Eq. 12 can no longer be ignored when constructing the linearized dynamics of a 
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two-dimensional homogeneous condensate and this leads to a contribution which breaks the 
symmetry of the two polarizations of the magnon modes describing the magnetization dy- 
namics of the condensate. Altogether the strong variation of post-quench dynamics with 
changes in the system geometry provides a compelling signature of dipole-dipole interactions 
that may be studied in future experiments. 

VII. NUMERICAL METHODS AND RESULTS 

Having set up our theoretical model we now turn to the numerical techniques and quan- 
titative results based on these ideas applied to a model spinor condensate with parameters 
motivated from experiment. As previously discussed, the calculation of correlation func- 
tions within the TWA requires the time evolution of an initial state which is comprised of 
an initial mean field state with random fluctuations added to it. The initial wavefunction 
of the rriz = condensate is determined by minimizing the total energy via conjugate gra- 
dient minimization assuming q oo. Time evolution according to the GPE is determined 
numerically by the 6-th order Runge Kutta methocl^^ with periodic boundary conditions 
in space. The kinetic energy is computed by Fourier transforming each component into 
momentum space. The dipole-dipole interaction kernel, K{r,r'), has the properties both of 
being long-ranged and also of being singular at short distances. Therefore it is necessary to 
regularize and truncate K{r, r') in real space before calculations are performed in Fourier 
space to avoid interaction between inter-supercell periodic images as discussed in Appendix 
B. In calculating K{r, r') for use in the solution of the full GPE, we neglect the variation of 
the condensate thickness {RTF,y) along the x direction. We have checked that this approxi- 
mation doesn't significantly affect our results when D is along the long axis of the trap i.e. 
z, as is the case in experiment. 

For the calculations reported we use a time step of 3.5 /is and a grid spacing of 0.5 /im 
and the results are found to be converged with respect to these parameters. In addition, 
the total energy of the system remains conserved to a certain error tolerance in the time 
evolution. It is also verified that the total magnetization along D is a conserved quantity in 
the absence of dipole-dipole interactions. 

The trap geometry for our calculations is taken to be similar to experiment- such that 
the TF radii of the condensate are -Rtf,x = 20 /im, Rtf.v = 1.6 /im. Given that the relevant 
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lengthscale for spin dynamics (27r/ ^j2M\(>^(ri^ ^ 2/im) is much smaller than the z length 
of the condensate {Rtf,z = 200 /xm), here we treat the system as unconfined along z, with 
periodic boundary conditions over a 90 fim length. The peak three-dimensional and two- 
dimensional densities are taken to be n^u = 2.5 x 10^^/cm^ respectively. The strength 
of the spin-dependent part of the contact interaction has been inferred previously from 
molecular spectroscopy2^i^ and from spin-mixing dynamics^^i^^. According to these works 
C2(r = 0) = 0.8 |c2|n3_D is predicted to lie between /i x 6 Hz and h x 8 Hz, corresponding to 
1.1 aB < Aa = (flo — 0,2) < 1.9 as where Aa is the difference between the s-wave scattering 
lengths for the spin-0 and spin-2 channels and is the Bohr radius. This variation in the 
rate of spin-amplification makes it difficult to estimate how close the estimated initial noise 
from experiment is to the quantum limit. 

We find our results to be in qualitative agreement with experiment and previous the- 
oretical calculations. In particular, we find that the average magnitude of the transverse 
magnetization grows exponentially from a small value to a much larger value (Fig. 1) with a 
time-constant that is relatively insensitive to the quadratic Zeeman shift q. The calculated 
domain structure and magnetization correlations match with those observed experimentally, 
and in previous calculations-, the characteristic domain size increasing with q (Fig. 4). How- 
ever as seen from Fig. 5 our calculations somewhat underestimate the domain size for the 
larger of the experimentally measured values of Aa. This discrepancy between theory and 
experiment is reduced on using the smaller of the measured values of Aa. Thus the differ- 
ence between theory and experiment could be the result of an error in the experimentally 
measured value of the spin-dependent contact interaction or quantum and thermal effects 
of interactions that are not contained in the Gross-Pitaevskii equations. The introduction 
of a dipole-dipole interaction introduces a weak dependence of the average local transverse 
magnetization G{0) on the quadratic Zeeman shift. 

Despite the qualitative agreement between the homogeneous 2D condensate calculation- 
and the current results, we find quantitative differences between the results of the homoge- 
neous case and the calculations including the trap and dipole interactions that are important 
for comparison to experiment. We discuss several of these differences below. 
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FIG. 4: Correlation function G{5r = zz) at t = 100 ms for a spinor condensate calculated using the 
method and geometry described in Section VII for various quadratic Zeeman shifts. The correlation 
function plotted along the length of the condensate shows decreasing domain size with decreasing 
quadratic Zeeman shift. 

A. Effect of the trapping potential 

Similar to previous theoretical work'', we find a significant effect of the inclusion of the trap 
on the spin dynamics in the parameter regime corresponding to experiment. The external 
trapping potential along the width of the condensate, which is accounted for in our numerical 
calculations, is found to slow the growth of the transverse magnetization in the condensate 
significantly as seen in Fig. 7, when compared to the quasi-two-dimensional homogeneous 
case without a trapping potential along the width. This slowing down can be understood 
intuitively from the fact that the trap causes the density away from the center of the trap 
to be lower than at the center of the trap. Consistent with previous theoretical work^, the 
density reduction away from the center of the trap also affects the spatial structure of the 
correlations observed and the trap is found to suppress the formation of structure in the 
radial direction, as seen in Fig. 5. 
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FIG. 5: Planar correlation function G{5r) at t = 87 ms for a spinor condensate calculated using 
the method and geometry described in Section VII for q/h = 2 Hz. The one dimensional plots 
shown in the lower half represent sections of the two dimensional plots above through the center of 
the condensate at x = 0. The two dimensional plots are 100 along the z direction and 40 
along the x direction. The above plots compare the calculated correlation structure on the left to 
the experimentally measured structure on the right^. 

B. Effect of the dipole interaction 

As discussed in Section VI, dipole-dipole interactions reduce the rate of domain formation 
in the case where the magnetic field is aligned along the z-direction, which is the long axis 
of the condensate. However, as seen in Fig. 5, for the parameters of the calculation, which 
are taken to be the ones relevant to experiment, the effect of the dipole-dipole interaction 
on the average transverse magnetization G{0) turns out to be small because the Fourier 
transform of the dipole interaction kernel K almost vanishes at the lengthscale of domain 
formation. That is, the spin healing length being nearly equal to the narrower condensate 
thickness, the dominant length scale for domain formation coincides with the cross-over 
between the 2D and 3D forms of the dipole interaction. Despite having a negligible effect 
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FIG. 6: Spatial structure of the complex transverse magnetization F±(r) in the x — z plane at 
t = 87 ms for q/h = 2 Hz. The upper panel of figures consist of random instances of experimentally 
measured spin textures^ while the lower figures are random instances of calculated spin textures. 
Each figure in the panel is 100 /xm long and 40 //m wide. The complex number F±_ is represented 
as shown in the color wheel in the inset of the figure. Domains are seen to have limited structure 
along the width of the trap and localized at the center of the condensate. 

on domain formation in the longitudinal direction, dipole-dipole interactions are found to 
suppress domain formation along the radial direction. 

Yet as discussed in Section VI, other experimental geometries, i.e. orientations of D and P 
away from the z axis, are expected to show more prominent dipolar effects in the spontaneous 
formation of magnetization. We explored this possibility numerically. The magnitude of the 
magnetization variance G(0) indicated by such calculations is shown in Fig. 6. One can see 
that the rate of growth of transverse magnetization is significantly enhanced with D and P 
pointing along the y-direction compared to other orientations. 

VIII. CONCLUSION 

In our analysis we have studied a realistic quantum Hamiltonian model for a quasi-two di- 
mensional spinor condensate, including the effects of the trap and dipole-dipole interactions, 
in order to make a quantitatively accurate prediction of the contribution of intrinsic fluctua- 
tions to the symmetry-breaking domain formation. Similar to previous studies, the inclusion 
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FIG. 7: Dependence of the evolution of the variance of the transverse magnetization for the dipole- 
precession axis D ahgned along x,y and z directions for q/h = 2.0 Hz. The growth rate of transverse 
magnetization is found to be enhanced significantly for D aligned along y as compared to D aligned 
along z which yields results close to the case without dipole interaction. 

of the trapping potential was found to reduce the rate of structure formation because of a 
reduction of the average density. The dipole-dipole interaction, which is known to have a 
dominant effect on the long term structure formation in spinor BECs, was found to add an 
effective non-local contribution to the spin-dependent part of the interaction in the spinor 
condensate. The non-local nature of the spin-spin interactions couples the spin structure 
formation dynamics to the direction of the spin polarization. Even though dipole-dipole 
interactions are found to affect the spin dynamics weakly when D and P are polarized along 
rriz = 0, we find that dipole interactions significantly enhance the rate of domain formation 
when these vectors are polarized along rriy = 0. Moreover, dipole-dipole interactions were 
found to split the degeneracy of the two polarizations of the magnon modes in the case 
where D was orthogonal to P. This spin-polarization dependence of the domain formation 
rate leads to a direct way to observe experimentally the role of dipole-dipole interaction on 
spinor dynamics. 

A quantitative prediction of the magnitude of the structure formation also requires the 
inclusion of effects from non-linear interaction terms and thermal effects. On a preliminary 
examination one would have expected thermal effects to be significant since the kinetic 
temperature of the condensate is much larger than the spin-mixing energy scale. However 
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FIG. 8: Summary of effects of various factors on the evolution of the magnitude of transverse 
magnetization. Both dipole-dipole interactions and the external trapping interaction are found to 
reduce the growth rate of the transverse magnetization. 

we found that the coupling of phonons to spin fluctuations is small, leading to a separation 
of the low temperature spin dynamics from the high temperature phonon dynamics. In 
Section IV, we studied the effects of the non-linear interactions within the standard TWA 
and found that non-linearity effects lead to saturation of the transverse magnetization at 
long times. We expect the TWA to be a reasonably accurate description of the spin-spin 
correlations of a spinor BEC since it was found to yield results in good agreement with exact 
diagonalization calculations our spinor BEC model in the single mode regime. 

Despite our effort to include the effects of the trapping potential, dipole-dipole inter- 
actions, non-linearities and finite temperature to develop a quantitative understanding of 
the magnitude of domain formation, we found in Section VII that the uncertainty in the 
magnitude of the spin-dependent part of the contact interaction prevents us from making 
a quantitative comparison of the magnitude of the domain formation with experimental 
results. Such a quantitative comparison between theory and experiment is critical for the 
determination of the contribution of intrinsic quantum fluctuations to the domain formation. 
One possible experimental approach to resolving this problem is to determine in a direct 
way the gain of the spinor BEC in the experimental geometry by studying the dynamics of 
the magnetization of the condensate following an initial microwave pulse. Such experiments 
in conjunction with quantitative calculations might make it possible to determine better the 
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importance of intrinsic quantum fluctuations to symmetry breaking dynamics. 
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APPENDIX A: EIGENMODES FOR THE BOGOLIUBOV TRANSFORMA- 
TION OF INHOMOGENEOUS DIPOLAR CONDENSATES. 

Here we give explicit expressions for the eigenmodes and eigenfrequencies of a dipolar 
ferromagnetic spinor BEC for positive quadratic Zeeman shifts. As discussed in Section 
VI, the inclusion of dipole-dipole interactions requires the generalization of the local spin- 
dependent coupling constant to C2*^(ri, r2). The eignmodes and eigenfrequencies that we 
define are strictly vahd when the Hermitean Hamiltonian Hq — —■^\/'^ + q{t)+fj, + con{r) + 
Vtrap is positive definite. In this case the eigenmodes of the condensate are given by 



where E'^ and ^„ are defined to be eigenvectors and eigenvalues of the Hermitean operator 
Hq^^{Ho + 2cf^)HQ^^ and a — ±1. The dual modes then follow to have the form 



In the case of negative quadratic Zeeman shifts, such a Hermitean eigenproblem cannot 
be constructed since some of the frequencies in this case are neither purely real nor imag- 
inary. This can be verified by introducing a weak periodic potential to the homogeneous 
ferromagnetic Bose gas and diagonalizing the problem using degenerate perturbation theory. 

APPENDIX B: REGULARIZING THE DIPOLE POTENTIAL 

In order to perform numerical calculations, even semi-analytic calculations where the 
Fourier transform of the three-dimensional dipole interaction kernel, K(R — R') = 



(Al) 
(A2) 



(A3) 

(A4) 
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|R~^R''[^~"^ ' needed, one needs the integral involved in the Fourier transform 
of the kernel to be well defined. The full 3D Fourier transform of the dipole interaction 
kernel may be calculated analytically, but to obtain converged results for the spin dynamics 
it is necessary to truncate the long-ranged dipole interaction between periodic images of the 
system which emerge when using Fourier techniques to do such calculations. For numeri- 
cal convenience we imagine that the dipole density can be expanded in terms of a possibly 
overcomplete set of functions i.e 



(grid) ^ 



E„p(R 



(grid) ^ 



(Bl) 



Such an expansion allows us to represent a function 0(R) which is smooth on the scale of the 



width of p(R) by its value 0(Rn^"'^^) on a discrete grid of points R^f"'''' with a grid spacing 
that is smaller than the width of p(R). Integrals of the kernel -R'(R) of interest are given by 



{grid) 



J 0(Ri)0(Ri)K(Ri - R2) = ^ 



j dR[ dR2p(R; - Ri)p(R2 - R2)i^(Ri - R2)- 



Taking the smoothing function to be p(R) 
given by 



rfRic;R2p(Ri - R)p(R2)i^(Ri - R2 
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the averaged kernel g is 
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(B5) 

The regularized expression, unlike the original dipole interaction kernel -ft'(R), vanishes for 
small R and approaches the regular expression for large R as expected. 

Using the regularized kernel g(R), the dipole interaction can be calculated on a real- 
space grid with grid spacings smaller than the width of the smoothing profile p(R), in a 
way so as to avoid interaction between periodic images. The Fourier transform for the two 
dimensional kernel with a Gaussian profile given Eq. [15] in Section VI can be derived by 
applying a Fourier transform to g(R) restricted to the 2D plane. 
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